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I. INTRODUCTION 

In the past few years, numerical techniques in fluid dynamics have matured to 
the point that meaningful simulations of three-dimensional transitional or turbulent 
shear flows bounded by a wall are possible. At the same time, availability of 
advanced computers, ovef and above a CDC 7600-class machine, has allowed roughly a 
doubling of the spatial resolution in each direction, greatly enhancing the quality 
of the results. We are in a position now to expand our set of objectives for flow 
simulations in the next decade to the following: 

1 . Increase our understanding of the physics of turbulence 

Provide statistical and other information useful to turbulence modelers 

3. Provide a testbed for the development of improved numerical methods and 
guidance for improved modeling of fine-scale effects. 

Objectives (1) and (2) are, of course, in common with the goals of most experi- 
mental (and purely theoretical) investigations of turbulent flows. The numerical 
simulation, however, provides a different way to look at a turbulent flow; the results 
of computer simulations are complementary to those of laboratory experiments. In the 
computer experiment, the corresponding laboratory time of the simulation is relatively 
small, but complete flow-field information is available throughout that time interval. 
The third objective is related to the desire to improve future simulations and to 
make it possible to Investigate more complex flows. Advanced computers, as they 
become available, will also help in this regard. 

In this paper we report on three recent simulations, using the ILLIAC IV com- 
puter, that begin to meet the above objectives: 

1. Study of the evolution of a wave-like disturbance in a laminar boundary layer 
corresponding to vib rating-ribbon experiments. 

2. Study of the evolution of a spot-like disturbance in a laminar boundary layer. 

3. Investigation of turbulent channel flow. 

In particular, we will concentrate on the turbulent structures observed in these 
flow simulations and comment on cheir role (1) in the dynamical process, for example, 
the growth and maintenance of turbulence, and (2) as contributors to turbulent statis- 
tics of interest to turbulence modelers. As will be shown, these persistent flow 
structures include streamwise and vertical vorticity distributions near the wall, 
low-speed and high-speed streaks, and local regions of intense vertical velocity. 
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The problem of representing the large range of physical scales in a computer simula- 
tion will also be discussed, along with some prospects for future improvement. 

Quite different numerical techniques are employed in the three simulations. The 
Urst and third simulations listed above use an Eulerian grid with pressure and veloc- 
ity as dependent variables and periodic boundary conditions in two directions, whereas 
the second simulation uses a Lagrangian vortex method in a semi-infinite domain. The 
t irst two simulations studied the evolution in time of an initial disturbance and 
b.i.-e no explicit model to account for "subgrid-scale" turbulence — fine-scale 
structures not resolvable by the computational mesh. Consequently, in these simula- 
tions the evolution time during which meaningful results were obtained was limited, 
due to the importance of subgrid effects at later times. In the study of turbulent 
channel flow, however, the objective was to simulate a fully developed turbulent flow 
having a range of scales much larger than the computational mesh would allow. Thus, 
the use of a subgrid model was mandatory. In Section II, further details of the 
simulations are given along with major findings regarding the observed turbulent 
structures. The problem of accounting for turbulence scales not resolved or repre- 
sented by the computational mesh — subgrid or supergrid modeling or both — is dis- 
cussed in Section III, and some possibilities for improvement are suggested. In 
Section IV, a summary of the results and conclusions are presented. 


II. RESULTS OF THREE NUMERICAL SIMULATIONS 


A. Flat-Plate Transition 

In a series of computer experiments Wray and Hussaini [1] (hereinafter referred 
to as WH) studied the instability of a flat-plate boundary layer in three dimensions 
by solving numerically the incompressible Navier-Stokes equations. In another three- 
dimensional numerical study, Orszag and Kells [2] investigated transition of plane 
Poiseuille and plane Couette flow. The ILLIAC computer code of WH is capable of 
simulating a wide range of boundary-layer stability problems. In this particular 
investigation they simulated the vibrating-ribbon experiments of Kovasznay et al. [3] 
by a careful choice of initial conditions for the three-dimensional perturbation to 
the Blasius profile. Two oblique (three-dimensional) eigensolutions of the Orr- 
Sommerfeld equation were added to a two-dimensional eigensolution to form the initial 
perturbation. To correspond with the laboratory experiment the streamwise and span- 
wise wavelengths of the perturbation waves were 25.2 6* and 30 6*, respectively, where 
6* is the initial displacement thickness of the laminar boundary layer. A schematic 
of the numerical simulation of WH is shown in Fig. 1. Wray and Hussaini used Fourier 
expansions (spectral method) in the streamwise (x) and spanwise (z) directions with 
64 points each and second-order finite differences in the normal (y) direction with 
128 points on a stretched mesh. Fourier transformation of the momentum and mass 
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FIG. 1 Simulation of a wave-like disturbance in a laminar boundary layer. 


conservation equations in x and z and some manipulations to minimize storage 
requirements yield the following system of equations for every (k^.k,): 
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Equations (I)~(4) give the time evolution of the Fourier transformed velocity compo- 
nents u, v, w and a pseudopressure q = p + v 2 . Here k 2 = l^ 2 + k z 2 . In comput- 
ing the nonlinear terms above, WH zero any aliasing contributions by the 2/3 rule, 
yielding a calculation that is equivalent to a Galerkin spectral method with 22 wave 
numbers in the x and z directions. 

The use of periodic boundary conditions in the streamwise direction leads to an 
evolution of the disturbance in time rather than space and to a slight departure from 
the Blasius profile. However, this assumption allows WH to apply 64 grid points 
(22 wave numbers) to each wavelength of the disturbance, thereby achieving high reso- 
lution of streamwise variations in the flow field. In the laboratory experiment, the 
evolution of the disturbance occurs in space over a distance of about 10 wavelengths. 
Thus, in the laboratory a particular structure at one streamwise location is a 
neighbor to the same structure located 1 wavelength downstream, but in a slightly 
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advanced state of evolution. Similarly, an earlier state of that structure Is situ- 
ated I wavelength upstream. As found by WH, however, the difference between the spa- 
t tul and temporal developments in the linear part of the evolution (exponential 
growth) is quite small as measured by the amplitudes and phases of the eigenmodos it' 
the two problems. Apparently the nonlinear evolution is very similar as well, -judgin} 
b ; the excellent agreement between the two experiments. In Fig. 2 are shown contour: 
•i ni/Dy for both the laboratory experiment and the computer simulation at the one- 
spike stage — the spike being the kink in the high amplitude 8u/3y contours [1], 
i'he contours are taken in the x-y plane containing the maximum initial perturba- 
tion the peak plane. At this point the time of evolution in the simulation is 
340 5*/U w . The agreement is quite good as is the comparison of the corresponding 
normal velocity contours given in Fig. 3. The maximum vertical velocity is roughly 
0.07 U M at this stage. The boundary-layer thickness, 6, is approximately 3 times 
the displacement thickness. 




FIG. 2 Contours of 8u/9y in peak plane — one-spike stage: (a) experiment [3]; 

(b) numerical simulation of WH. 
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FIG. 3 Vertical velocity in peak plane — one-spike stage: (a) experiment [3]j 

(b) numerical simulation of WH. 

As nonlinear effects become increasingly significant, the evolution of the com- 
putational disturbance reaches the two-spike stage at t = 360 6*/Uc<,. Contours of 
9u/3y and normal velocity at this stage are shown in Figs. 4 and 5, respectively. 

Note that again there is very good agreement on the structure of the disturbance in 
the peak plane. The relative shift of the structures along the abscissa is irrele- 
vant [1], but other discrepancies are evident. For example, the computation is pro- 
ducing peak values in v ( = 0.2I U^) that are higher than the experiment. Considering 
the large streamwise gradient of v in the vicinity of this maximum, it is quite 
possible that this type of discrepancy is due to lack of resolution or to ensemble 
averaging or both in the laboratory experiment. On the other hand, the appearance of 
oscillatory fine-scale features in Fig. 4 indicates that the computation may be start- 
ing to have a resolution problem also. A three-spike stage occurs in the simulation 
at t = 365 6*/U 00 . Soon after this time significant errors begin to appear due to 
the fact that the energy of the disturbance is trapped to remain within the computa- 
tional wave numbers. 
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FIG. 4 Contours of 3u/3y in peak plane — two-spike stage: (a) experiment [3]; 

(b) numerical simulation of WH. 

Thus, in a remarkable computational effort WH have demonstrated the capability 
of simulating a complex, three-dimensional transitional flow in an "honest" calcula- 
tion — no subgrid modeling was used. In the particular case that was stimulated the 
results appear to be valid up to and including the two-spike stage when v max = 0.21 U w . 
After that stage, presumably, scales with wavelengths <306*/22 in the spanwlse 
direction or <256*/22 in the streamwise direction become important. The vertical 
resolution appears to be quite adequate throughout the simulation. 

B. Spot-Like Disturbance in a Laminar Boundary Layer 

In another numerical Investigation of boundary- layer mechanics, Leonard [4] 
studied the growth of localized disturbances in a laminar boundary layer as they 
evolve into spot-like turbulent structures. As in the study of WH, these are computer 
experiments, but using a quite different numerical method. Vortex filaments are used 
to represent the three-dimensional, perturbed vorticity field within the spot. Each 
filament is defined by a three-dimensional space curve x(£J,t), circulation r, and a 
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FIG. 5 Vertical velocity in peak plane — two-spike stage: (a) experiment [3]; 

(b) numerical simulation of WH. 

core size, a, which parameterizes an assumed Gaussian vorticity within the core. 

Thus, the perturbed vorticity field is given by a summation over filaments, 

Ji'(x*t) - f exp[-|x - x 1 | 2 /0 2 ]/(/Fa) 3 (3x i /3?)d5 . (5) 

The vortex filaments represent a continuum of vorticity which initially is in 
the form of an infinite sheet of finite thickness, that is, a boundary layer* As 
indicated in the schematic (Fig. 6) the disturbance in the simulation is created by a 
local distortion of the vortex filaments at t ■ 0. For numerical purposes each tube 
of vorticity in the boundary layer is decomposed into its straight, unperturbed con- 
tribution and into a loop representing its contribution to the perturbed vorticity 
field, as shown in Fig. 7. The velocity field is decomposed, accordingly, as 

u(x) - U(y)e x + u' (x) , (6) 
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FLOW TANQENCY AT V " 0 
u(k) -» U(y) # x |k| •*“ 


where U(y) Is the laminar profile 
and ju* is given by a sum of Biot- 
Savart Integrations over all the 
filament curves, plus their 
Images. The curves are marked with 
a sequence of node points that are 
tracked In a Lagrangian reference 
frame, according to the following 
equation of motion: 



(Xj_~ x^xQXj/H^dS 1 
[lit “ x j 1 2 + 01 C *^ 2 + 0j z )] 3 /c 

+ e x U(y) + image contributions. 

(7) 

The equation is seen to be the Blot- 
Savart induction law for an infini- 
tesimal filament, but corrected in 
the near-field for finite core size. 
Using a = 0.2065 we recover the 
correct speed of a ring vortex in the 
limit a « ring radius. Up to 
numerical truncation error, the 
method simulates the inviscld dynam- 
ics of the vortlcity field and there- 
fore satisfies the Kelvin and 
Helmholtz theorems. For a more complete description of the method including simula- 
tion of viscous effects see [5]. 

Blot-Savart integration yields exact boundary conditions at infinity in all 
directions. The image contributions ensure tangency of the flow at the wall. In this 
initial study we ignore the generation of vortlcity at the wall boundary due to the 
disturbance, that is, we use the inviscid boundary condition at the wall for u' . 

This appears justified because the relatively high Initial disturbance level 
(0.07 - 0.10 U„) and the initially small length scale of the perturbation “0(5*) 
(corresponding to laboratory experiments on turbulent spots [6,7]) leads to the rapid 
evolution of a sizable, strong disturbance at the end of the simulation, t ~ 506*/U ro . 
The viscous length scale for this time, assuming Re^ = 1000, is only 0.36*. 

Computational filaments are added automatically as the disturbance grows in size. 
At the end of the simulation, about 8000 points represent one half of the spot 


FIG. 6 Simulation of a spot-like distur- 
bance in a laminar boundary layer. 



FIG. 7 Decomposition of vortlcity field. 
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(symmetry of the flow about the mldplana in enforced) , The computation could contlnuo 
beyond this time, but (1) some of the filaments are quite stretched at this point, 
indicating that remeshing is required, and (2) the cpu time per time step is quite 
large at this stage being proportional to the square of the number of computation;*! 
points for vortex methods. Total run times are =2 hr on the ILklAC. 

In Fig. 8, the space curves defining the vortex filaments originally in a single 
layer (out of a total of six representing the boundary layer) a in shown at three dif- 
ferent times. The spreading of the disturbance and the increasing deformations suf- 
fered by the vortex lines are quite evident. In Fig. 9, vertical and streamwise vor- 

ticity contour plots at a late time are shown in horizontal (x-z) and spanwise- 

vertical (z-y) planes, respectively. Note the development of Intense streamwise and 
vertical vorticity generated in the flow in agreement with the pattern of vortex llneb 
displayed in Fig. 8, The structures are cell-like, 0(6*) in the spanwlse and vortical 
directions, but elongated in the streamwise direction. For example, the region of 
negative w y just above the x-axis in Fig. 9a together with its antisymmetric coun- 
terpart (not shown) just below the axi" wouU constitute one cell approximately 46* 

in spanwise extent. Velocity contour plots are shown in the same (x-z) and (y-z) 
planes in Figs. 10 and 11, respectively. The presence of streaky structures in the 
horizontal plane is compatible with experimental observations of passive particle 
motion in turbulent spots [7] and turbulent boundary layers [8). The maximum normal 
velocity is quite high at *0.4 U w . 

Most laboratory experiments on turbulent spots report ensemble-averaged measure- 
ments several hundred 6* downstream of the Initial disturbance. A typical stream- 
wise velocity history plot from Cantwell et al. {7] measured in the midplane of the 
spot is shown in Fig. 12a. Note that compared with the laminar profile there Is a 
velocity excess near the wall and a defect away from the wall — typical of a turbulent 
boundary-layer profile. The experiment of Wygnanskl et al. [6] produced similar 
results. However, as indicated in Figs. 10 and 11, the particular numerical simula- 
tion presented above yields a velocity excess in the midplane for all heights above 
the boundary. On the other hand, if the velocity field from the numerical simulation 
is span-averaged over one cell of width 46* centered on the midplane we obtain the 
velocity plots shown in Fig. 12b. (TO obtain better statistics, the results shown in 
Fig. 12b are based on the velocity field at a single time late In the simulation. 

Thus, decreasing x is used along the simulation rather than increasing time.) Note 
the close similarity in the behavior of these computational results averaged over one 
cell in the span direction at t - 43 6*/U M to the ensemble-averaged laboratory mea- 
surements of Fig. 12a, taken at a fixed streamwise location, at later times 
(t > 1500 6*/U 00 ) . The good comparison leads to the speculation that the internal 
structure of a single realization of a large turbulent spot (>1006* in streamwise 
extent) is dominated by the presence of many cell-like substructures such as those 
revealed in the numerical simulation. (For experimental support of this speculation 
see I. Wygnanskl, these proceedings.) 
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FIG. 9 Contours of vorticity for 
y/6* - 0.66; (b) a) 


lCt/6* >= 41; (a) 


0), 


. , „ - ■ J y in a horizontal plane, 

>x in a vertical (z-y) plane x/<5* « y 25. , positive contours- 

negative contours; contour levels +0.05, +0.15, +0.25 . . .of 8u/3yjwall. ’ 


Three-dimensional vortex simulations have therefore demonstrated the rather 
rapid, nonlinear evolution of small localized disturbances in a laminar boundary layer. 
The boundary layer appears to deform into a corrugated surface with significant amounts 
of streamwise and vertical vorticity - leading to strong negative u'v correlations 
and, hence, intense turbulent mixing. Velocity perturbations induced in the nee-by, 
relatively undisturbed vorticity layer lead to sufficient deformations in the vortex 
lines to produce local nonlinear growth. Thus, the spot appears to grow in the span- 
wise direction by "transverse contamination," in agreement with recent experiments 
using dye interjection [9], which demonstrated that turbulent convection is not 
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FIG. 10 Velocity contours in the horizontal plane of Fig. 9a. 

Contour levels ±0.02, ±0.06, ..0.10 . . . of U^. 

responsible for the spanwise growth of a spot. Further study is needed to simulate 
properly larger turbulent spots or a collection of spots and the decay of turbulence 
rearward of the spot. 

C. Turbulent Channel Flow 

In our third example of a wall-bounded shear flow we discuss the numerical simu- 
lation of turbulent channel flow by Kim and Moin [10] (hereinafter referred to as KM). 
As shown in Fig. 13, a 64 3 grid is used to simulate flow at Re = 13,800 with periodic 
boundary conditions in the spanwise and streamwise directions and the no-slip condi- 
tion at the channel walls. In an earlier investigation of turbulent channel flow, 

Moin et al. [11] used a 16x64x16 grid with the same boundary conditions. Because the 
laboratory flow at this Reynolds number contains a significant amount of turbulent 
energy at length scales smaller than the computational grid, KM used the filtered 
Navier-Stokes equations to determine the dynamics of the turbulence field that is 
resolvable on their mesh (large eddy simulation technique). Thus, a flow variable f 
may be decomposed as follows: 

f - f + f ' , (8) 
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where f is the large-scale field, 
defined as the result of a filtering 
or smoothing operation of f, and f' 
is the subgrid field. 

After applying the filtering 
operation to the incompressible 
momentum and continuity equations, 
the governing equations for the fil- 
tered field may be written 


3u 

Tt 

and 


+ uxu 


dP 


-Vp* + ^ e x - V*jr + vV ? u 

(9) 


V*u = 0 , 


00 ) 


where dP/dx is the constant mean 
pressure gradient, iu is the large- 
scale vorticity field, and the 
subgrid-scale stress tensor, jr, 
(having zero trace) has components 


defined by 



T as 

ij 

8 « - W 3 

(11) 

and 



“ij " u i u j 

+ ulu. + u.uj . 
j i 3 i 

(12) 


FIG. 12 Streamwise velocity: (a) ensemble- 

average velocity records [7], x/6„* “ 1600; 
(b) cell-averaged streamwise velocity from 
numerical simulation of [4]. 



U c H/2 

Re-- 2 — =13800 
GRID:*’ 64 X 64 X 64 
BOUNDARY CONDITIONS: 
u — v — w — 0. AT y — i H/2 
PERIODIC IN x AND z 


The quantity p* is given by 

P* - P/P + 1/2 + R^/3 . (13) 


The system of equations (9) and (10) 
becomes deterministic once x_ is 
specified. Kim and Moin use an eddy 
viscosity model, 

T ii * - 2V T S iJ ' < W > 

where 


3 ij 


" 2 ^x. 9x i y * 


(15) 


and 
model, 


is given by the Smagorinsky 


<v>V 


S..S,. . (16) 

ij ij 


FIG. 13 Numerical simulation of turbulent 
channel flow. 
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The constant Cg is taken to be O.L and A is a measure of the grid resolution and 
is taken to be 

A ■ (Ax Ay Az) l / 3 [1 - exp(-y + /50) ] (17) 

where the exponential function damps the eddy viscosity near the wall. 

The channel height is H. The lengths -vjf computational box in the streamwise 

and spanwise directions were chosen to be L x « a nd L z * 2ttH/ 3, respectively, 

based on two-point correlation measurements in laboratory experiments and the desire 

to capture all the important large eddies in the flow. The nondimensional grid spac- 

ings that result are Ax - 63 and Az “42, expressed in wall units. (The channel 

height H is 1280 wall units.) The grid spacing In the normal direction varies from 
4* 4 

Ay z 2 near the wall to Ay 2 50 at the middle of the channel. The need for sub- 
grid modeling is clear. Even near the wall, where grid is finely meshed in the normal 

direction, the horizontal spacing is not sufficient to resolve the experimentally 
observed wall streaks with a mean spanwise spacing of \ z ~ 100. Kim and Moin find, 

however, that wall streaks do appear in their simulation, but at a larger scale, the 

streaks apparently accommodating themselves to the available resolution. 

The agreement of KM's simulation 
with laboratory measurements is 
remarkably good. After a time of 
evolution of 15.4H/U C (U c is the 
mean velocity of the center of the 
channel) KM obtain the mean-velocity 
profile (horizontally averaged) 
shown, along with experimental data, 
in Fig. 14. Near the wall the simu- 
lated profile is somewhat low; away 
from the wall a distinct logarithmic 
region is revealed with the correct 
slope (Kdrman constant) . The com- 
puted large-scale turbulent intensity 

profiles and the turbulent shear-stress profile are shown in Figs. 15 and 16, 
respectively. The ^u" 2 ^ 1 ^ 2 and (w 2 ^ 1 / 2 profiles are in excellent agreement with 
experiment at all distances from the wall. (Superscript double prime denotes devia- 
tion from horizontal average.) In the wall region (y + < 90) the resolved vertical 
turbulence intensity is lower than experiment, possibly due to the fact that a sig- 
nificant amount of ^ V 2 ^ 1 / 2 resides in the subgrid scales or that the subgrid model 
has deficiencies near the wall, as suggested by the deficit of the mean profile in 
this region. 

Contour plots of instantaneous Reynolds stress u"v reveal interesting phenom- 
ena. In a horizontal plane near the wall (y + ■ 16, Fig. 17) a number of hotspots 
occur where u"v is large and negative. Kim and Moin determined that a great 



FIG. 14 Mean-velocity profiles; — ©-, 
numerical simulation of [10]; experiment, 

, Re = 13800; , Re = 57000; , 

Re = 12200-61600. 
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FIG. 15 Comparison of the horizontally averaged resolvable turbulence intensities 
normalized by u-t with experimental data; ®, numerical simulation of [10]; experi- 
ment, , Re = 13800; , Re = 15200; , Re = 25000; , Re = 120000. 

majority of these active regions is associated with u" > 0, v < 0, or "sweeps." 

Away from the wall (y + = 90, Fig. 18) the hotspots are associated with "bursts" or 
u" < 0, v > 0. All these findings are in agreement with experiment [12]. 

Kim and Moin also investigated pressure-velocity gradient correlations in the 
region of y < 100. These correlations provide for energy transfer between the 
individual components and turbulent energy and, as such, are of great interest to 
turbulence modelers. The computed correlations are shown in Fig. 19. Their behavior 
for y + > 20 corresponds to the expected transfer from {u" 2 ) to (v 2 ) and (w 2 ). 

Close to the wall, however, KM show a substantial loss of {v ? ) and a corresponding 
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gain of \w 2 ) through the action of 
these correlations. By studying the 
fluid motion near the wall using 
"snapshots" or instantaneous contour 
plots, KM determined that the energy 
transfer was due to high-speed par- 
cels of fluid approaching the wall 
and then spreading laterally ("splat- 
ting" effect). 

D. Comparison of the Simulations 


In Table I we compare a number 
of features of the three simulations. 
The natural length and time scales 
for the first two studies or the 
evolution of a disturbance in a 
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FIG. 16 Vertical profile of horizontally 
averaged resolvable turbulent shear stres 
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FIG. 17 Contours of u"v in the x-z plane at y 














TABLE 

I 



SIMULATION 

PARAMETERS 



Wave 

Spot 

Channel 


disturbance 

disturbance 

flow 

Grid 




NxXN y xN z 

64x128x64 

500 filaments 

8000 points (1/2 spot) 

64x64x64 

Mesh spacing 




Ax/6* 

0.39 

0.41 

— 

Ay/6* 

0 . 0078 *+■ °° 

0.41 

— — 

Az/6* 

0.47 

0.66 

— 

(wall units) 
Ax + 

I6 a 

I6 a 

63 

Ay + 

0. 31 *»• " 

16 

1.8 + 50 

Az+ 

19 

26 

42 

Turbulence volume 

0.3 x 10 9& 

0.4 x l0 9b 

12 x 10 9 

(wall units 3 ) 




Time step 




AtUoo/6* 

0.24 

0.49 

— 

At + 

0.38 a 

0.78 a 

0.61 

Number steps/run 

1500 

100 

2000 

ILLIAC time /run, hr 

12 

2 

12 


a Estimate using 6* + = (u^U^Re^ = 40; (6*/U m ) +s = (u x z /Uo, 2 )Re^ ~ 1.6. 
^Estimated vertical extent of turbulence = 250 wall units. 


laminar boundary layer are 6* and 6*/U„, respectively. In the simulation of turbu- 
lent channel flow, appropriate length and time scales for the near-wall region are 
v/u T and v/u x 2 , where u T is the wall friction velocity, u T = jA> 9 (y) 8 y| wa n< To 
compare the transition studies with the channel flow simulation, the scales of the 
first two studies were converted to approximate turbulent scales as follows. Assume 
that the laminar boundary layer completes its transition to a turbulent layer. Then 
the original laminar displacement thickness in wall units is given by 



1000. The time 


I 


% 


i 


f 

i 

f 

r 

\ 




6* n 


6*u t u T 


U* Re 6* 


s 40 , 

where we have used u T /U w * 0.04 for the turbulent layer and Rc^ » 
scale in wall units is found to be 



■> 1.6 . 


Note that the two transitional studies have comparable resolution in the horizontal 
directions, but the vertical grid of WH is much finer. Both are able to resolve 
structures having length scales 0(6*) in any direction, and,, as shown in earlier 
figures, many such significant structures were found. Recall, for example, the verti- 
cal velocity contour plot 
(Fig. 5) at the two-spike stage 
of the study by WH, exhibiting 
large peak values (5=0.21 U m ) 
and very large streamwise gra- 
dients. As shown in Fig. 20a, 
a similar structure with a peak 
v' of 0.42 U ro was found in 
the spot simulation. The 
results are shown in a z-plane 
(z « 0.8 6*) containing mostly 
upward fluid motion (bursting ?) 
similar to the peak-plane 
behavior observed by WH. Con- 
tours of the perturbed spanwise 
vorticity are shown in Fig. 20b 
in the same x-y plane at 
z = 0.8 6*. As might be 
expected, the region of large 
3v/3x is associated with a 
concentration of 

The time evolution of the 
maximum vertical velocities is 
shown in Fig. 21. Also shown is 
the range of peak values in the 




x/6* 


FIG. 20 Contours in an x-y plane from the 
numerical simulation of a spot-like disturbance: 
(a) vertical velocity; (b) perturbed spanwise 
vorticity. 
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channel flow simulation by KM In 
several y-z planes for y * 46. 
Through the final stages of evolution 
the peak vertical velocity in the 
transition study of WH is experienc- 
ing exponential growth at 15 time, 
the growth rate of the initial two- 
dimensional perturbation wave. 
Apparently this peak is associated 
with a small-scale (0(6*)) secondary 
instability of the larger (0(256*)) 
slower-growing primary Instability. 

The growth rate seen in the spot 
simulation is even larger and noL 
exponential, presumably due to non- 
linear interactions between small- 
scale perturbations. 

Note that in the channel-flow 
study a relatively coarse mesh in 
the horizontal plane was used, clearly 
dictated by the relatively large vol- 
ume of turbulence being simulated. 

A number of flow structures seen in 
the transitional studies that were 
(0(6*)) and smaller in the spanwlse 
and streamwise directions would not 
be resolved adequately by this mesh, 
(Recall that 6*^ ~ 40.) Nevertheless, KM see in their simulation considerable tur- 
bulent activity near the wall, as evidenced by Figs. 15-19. As was mentioned above, 
turbulent structures do appear (wall streaks, etc.) E but at a larger scale than in 
experiment, accommodating themselves to the available resolution. In Figs. 22 we see 
this effect illustrated in a contour plot of streamwise vorticity from KM's simula- 
tion. In Fig. 22a, about one seventh of the y-s computational plane is shown. Note 
that near the wall the spanwise extent of the cell-like structures is only a few times 

X 

the smallest available scale in the span direction (Az = 42). In Fig. 22b a blowup 
of one of the more active regions near the wall is shown. With one exception, the 
structure is similar to the one shown in Fig. 9b from the spot simulation. The verti- 
cal scales in Figs. 22a and 22b have been expanded by a factor of 4 over the horizon- 
tal scales so that the actual structure depicted in Fig. 22b is considerably flattened. 
For comparison, Fig. 23 shows Fig. q b from the spot simulation .sing the same scale in 
wall units of the channel-flow results. The smaller computational domain affords the 
luxury of seeing the fine-scale structures depicted in the figure. 
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FIG. 21 Maximum vertical velocities in a 
plane. 



FIG. 22 Contours of streamwise vorticity in a y-z plane — channel-flow simula- 
tion of KM: (a) one seventh of computational plane; (b) blowup of region indicated. 
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FIG. 23 Contours of streamwise vorticity in a y-z plane; spot simulation 

With scale in wall units. 


III. EFFECTS OF SUBGRID-SUPERGRID SCALE MOTIONS 


In the numerical simulation of fluid flows one can resolve only a limited range 
of length scales in the flow field. In most cases of Interest, turbulent flows will 
contain Important scales of motion outside this range,, A major problem in the simu- 
lation of turbulent flows, therefore, is to model the effect of these subgrid or 
supergrid components on the resolvable scales. The objective of this modeling effort 
should be to minimize the error on any desired statistical Information that is derived 
from the resolvable scale motions, including structural information. This problem 
appears particularly acute in the case of wall-bounded shear flows because of the 
importance of small-scale structures near the wall in the generation of turbulent 
energy and in turbulent mixing. In addition, the use of periodic boundary conditions 
will be a problem if scales larger than the computational domain are important. 

In the first two studies discussed in the previous section, which were concerned 
with the transition of a laminar boundary layer, there was no explicit attempt to 
account for subgrid-scale effects. In the study by WH of the evolution of a wave-like 
disturbance, a Fourier representation was used for the 'dependence of the flow variables 
in the streamwlse and spanwise directions. Thus, by the quadratic nonlinearity of the 
pressure and convective derivatives, the evolution of a given Fourier coefficient with 
two-dimensional wave number k 1 within the resolved range (k^n < | k x j < kfl, ax ) will 
contain contributions from other wave numbers kj and k 3 satisfying 

k 2 + k 3 = ki , 
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If nelthar of tha corresponding Fourier coefficients is negligible. If k., and k,j 
are in the resolved range there is no problem. If either or both are out of the 
resolved range an error is made by neglecting this contribution. The assumption of 
streatnwlse periodicity lead® to errors involving supergrid components. For the simu- 
lation of WH, this error (discussed lit Section XI) appears to be small. In flows 
where turbulent onergy is simply cascading to higher wave numbers near the cutoff 
^max ( e *8*> the decay of isotropic turbulence) the error involving subgrid components 
has the effect of turning off the cascade, resulting in a "pileup" of turbulent energy 
around |kj £ k,,,^ in the numerical simulation. Several schemes, based on a Fourier 
representation of the flow variable®, have been proposed to model the effects of the 
unresolved flow field (|k| > k,^) and, therefore, to allow the turbulent cascade 
process to pass energy on to the subgrid scales. See Leslie and Quarini [13) for a 
discussion of these possibilities. In wall-bounded shear flows, however, the produc- 
tion of turbulent energy may also occur near the cutoff wave number. This poses an 
additional challenge to the capabilities of a subgrid model. In this case the model 
must account for all the important unresolved Instability mechanisms that break up the 
wall shear layer and produce turbulence, 

Again, in the vortex simulation of a spot-like disturbance no explicit model for 
the subgrid scales has been used as yet. Fine-scales of turbulence are produced auto- 
matically in the simulation aB the space curves representing the vortex filaments take 
on highly convoluted configurations under the continual action of local velocity 
gradients. There is a lower limit to scale of turbulence that can be resolved, how- 
ever, From Hq. (5) we see that this scale is the core radius of the filaments, o, 

(Note that the Fourier transform of w contains the factor exp(-|k| 2 0 2 /4) so the 
maximum resolved wave number is "l/o.) In the vortex simulation discussed in this 
paper some of the computed space curves representing vortex filaments had become 
highly stretched and contorted by the end of the simulation. More computational node 
points would be required to accurately follow the subsequent motion of these filaments. 
fi -<*.bb costly, alternative way of proceeding with the simulation would be to remesh 
the vorticity field, perhaps by smoothing the filament curves and by changing the 
typology of some curves [5], Of course this procedure would be required repeatedly at 
later times as more fine scales reappear. The method appears to have promise as a 
subgrid modeling technique. Fine scales are allowed to develop in a natural way, 
retaining some interplay between the large and small scales. The method seems analo- 
gous to computing with a fine mesh in Eulerian calculation and periodically filtering 
out the high wave-number components of the velocity or vorticity field [14]. 

In the channel-flow simulation, KM use an eddy viscosity model for the subgrid- 
scale stresses. This model has proved adequate in less demanding situations, for 
example, the decay of homogeneous turbulence [15]. The model is successful in these 
cases in spite of the fact that it does a relatively poor job of estimating the actual 
point-wise subgrid-scale stresses computed by filtering the results of a fine-mesh 
simulation of a homogeneous turbulent flow [16,17]. Such a capability is, of course, 
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noc a necessary requirement for a subgrid model, but it would be more encouraging If 
tlm correlation between the actual and modeled stresses was higher than the W/, level 
determined by McMillan and Ferziger (17], A new model that does quite well in this 
respect has been proposed by Bardina et al, [18], The idea is to model the subgrid 
stresses directly in terms of the small-scale component of the reaolvable field. 

In the simulation of KM, there is every indication the the eddy-viscosity model 
It* quite satisfactory away from the wall (y + > 90). The mean velocity and the throe 
turbulence-intensity profiles agree well with experiment in that region. Near the 
wail, however, the mean velocity is low, as are the vertical velocity fluctuations, 
Apparently the eddy-viscosity modal, being essentially dissipative, is not able to 
represent or account for important turbulent production processes near the woll clrrt 
oecu< n the eubgrld scales or in the smallest resolved scales where the Interaction 
with the subgrid scales is of great importance. Evidence from the two studies con- 
cerning transition of a laminar boundary layer strongly suggests that that eventual 
breakup of the laminar layer involves structures that are 0(6*) in the shreamwisu 
and spanwise directions. In KM' s channel-flow simulation it is quite possible then 
that the lack of resolution in the x and z directions or other numerical effects or 
both inhibit the breakup of vortlcity layers at the wall until they grow excessively 
'hick. Indeed, in numerical experiments involving the filtering procedure, KM showed 
that substantial Increases in turbulent shear stress near the wall are possible for 
shore time intervals by modifying the treatment of the smallest resolved scales. It 
appears, therefore, that a subgrid model for the near-wall region should contain a 
destabilizing feature that would promote the breakup of vortlcity layers rather than, 
or in addition to, the dissipative aspect of the eddy-viscosity model. 

By using periodic boundary conditions in the horizontal directions, KM encoun- 
tered errors due to unresolved supergrid components of turbulence. The magnitude of 
this error is not known, but is estimated to be small based on experimental two-point 
correlations. 


IV. SUMMARY AND CONCLUSIONS 

Numerical simulations of wall-bounded turbulent shear flows have confirmed the 
existence, and illuminated some of the characteristics, of cell-like substructures 
near the wall containing substantial vertical and streatnwise vorticity. These struc- 
tures appear to be responsible for or intimately related to many features of turbulent 
flows observed experimentally — intense vertical mixing, low-speed and high-speed 
streaks, and, hence, the streaky structures seen in photographs of dye-particle 
contaminants . 

These structures were first revealed numerically in the simulation of turbulent 
channel flow. In this case, these structures seem to be related to strong (~u'v) 
hotspots and, ultimately, the mean-velocity profile since planar averages of u'v 
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yield Che turbulent shear-stress profile. In the case of a turbulent spot, span- 
averaging over one -ell produces a streamwiae profile qualitatively similar to 
ensemble-averaged measurements of experiment. 

In the simulation of the vibrating ribbon transitional flow all the features of 
the corresponding experiment were reproduced including the development of secondary 
instabilities with length scales 0(5*). At this writing, Wray is investigating the 
properties of a strong streamwise vortical structure that he has located off the peak 
plane In the simulation of WH. 

The impetus to imulate more complex wall-bounded shear flows in reasonable 
computer times means a continued sparseness of grid points for the near-wall region. 
More sophisticated subgrid-scale models will be required to account for the fine- 
scale structures discussed above that appear to be so important to the turbulent shear 
stress and turbulent mixing, in general, near the wall. Many interesting possibili- 
ties exist. We speculate that successful models will enhance the breakup of wall- 
bounded vorticity layers while perhaps retaining some aspects of an eddy-viscosity 
model. 
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